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We have studied by Molecular Dynamics computer simulations the dynamics of water confined 
in ionic surfactants phases, ranging from well ordered lamellar structures to micelles at low and 
high water loading, respectively. We have analysed in depth the main dynamical features in terms 
of mean squared displacements and intermediate scattering functions, and found clear evidences 
of sub-diffusive behaviour. We have identified water molecules lying at the charged interface with 
the hydrophobic confining matrix as the main responsible for this unusual feature, and provided a 
comprehensive picture for dynamics based on a very precise analysis of life times at the interface. We 
conclude by providing, for the first time to our knowledge, a unique framework for rationalising the 
existence of important dynamical heterogeneities in fluids absorbed in soft confining environments. 


I. INTRODUCTION 

Nanoconfinement and interfacial interactions pro¬ 
foundly alter thermodynamical, structural and transport 
properties of confined fluids. A molecular-level under¬ 
standing of these modifications is essential in a diver¬ 
sity of fields, ranging from biology [Q to nanoscience jj]- 
A situation of paramount interest is that of water 
adsorbed in natural or synthesized soft media, such 
as biomolecules, macromolecular gels, polymers, self- 
assembled and active materials. In particular, it is be¬ 
lieved that the function-structure relationship in these 
systems is primarily correlated to the amount of the ad¬ 
sorbed fluid and to its peculiar properties. Water in the 
vicinity of soft charged walls is expected to show sig¬ 
nificant deviations from the bulk behaviour, e.g., dis¬ 
tortions of the hydrogen bonds network and anomalous 
diffusion character, which might in turn profoundly im¬ 
pact the final morphology and functional behaviour of 
the confining media. A deeper insight of such complex 
fluid-matrix interplay can be obtained by using numeri¬ 
cal tools. Here, we use Molecular Dynamics simulations 
of a self-assembling model system to demonstrate and ra¬ 
tionalize the general existence of spatially heterogeneous 
dynamics in a fluid adsorbed in soft confining environ¬ 
ments. 

The generic situation of soft mobile confining media 
has been less explored than that of hard confining mate¬ 
rials. Detailed studies on confinement in highly stylised 
immobile matrices have been described, which aim at 
understanding the universal phenomena originating from 
both the detailed nature of the fluid/matrix interactions 
and the restricted character of the available space. These 
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include, together with many other systems, simulations 
of water in carbon nanotubes f3j, confined between two 
smooth hydrophobic plates [4] or patterned walls m , or 
in the case of patterned hydrophilic and hydrophobic 
walls [B]. In contrast, most studies on water confined 
in (disordered) soft matrices are very specific, and focus 
on particular confining macromolecules, without any con¬ 
cern for generality. These investigations include, among 
others, neutron scattering and molecular dynamics simu¬ 
lations of tau and MPB proteins ,7], dynamics of water in 
reverse surfactant micelles |8], or simulation of polymer 
electrolyte membranes [9, 10J. These molecules, be they 
of biological or synthetic origin, are very complex, and 
one finds extremely difficult to separate properties that 
are specific to the particular chemistry of the confining 
matrix, from those that are general and emerge imme¬ 
diately when water is in contact with any soft confining 
media. Still, in general we expect that the main fea¬ 
tures of water transport in these systems must conform 
to some universal description, dictated by the principles 
of statistical mechanics. 

It is well understood that water transport around pro¬ 
teins is highly anomalous EHH, with sub-diffusive wa¬ 
ter molecules playing a crucial role in the function of the 
protein. Notwithstanding this well-established observa¬ 
tion, the existence of sub-diffusion in synthetic materi¬ 
als is scarcely discussed. Surprisingly, most studies on 
synthetic materials focus on very complex simulations, 
with very detailed force fields, but they generally do not 
capitalise on this modelling accuracy to provide detailed 
and up-to-date analysis of the transport properties of the 
ab(d)sorbed fluid. For instance, the dynamics of water 
in surfactants is treated only in the reverse micelle region 
of the phase diagram HE and most studies consider wa¬ 
ter adsorbed in proton-conducting fuel cells ionomers as 
a Fickian diffusive system (9, Hfi] . Yet, the functional 
properties of these materials are also determined by the 
mobility of water confined in soft ionic nanochannels. 
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FIG. 1. Typical simulation snapshots of the ionic surfactants model studied here, at hydration levels A = 4,16,32, from 
left to right. Lamellar, cylindrical and micellar self-organised phases are evident. The hydrophobic section of the surfactant 
macromolecules is represented by the grey beads, the charged hydrophilic units by the green ones. Water molecules are indicated 
in blue, the hydronium cations in red. These configurations are examples drawn from the larger set of matrix templates we 
have considered in our investigation of the dynamics of water molecules confined in self-assembled soft environments. 


In addition, from the point of view of the elementary 
water-matrix interactions, ionomers, for instance, sub¬ 
stantially resemble biological matter: i) they both con¬ 
tain hydrophilic and hydrophobic parts, ii) neutral and 
charged groups, in) water is confined at the nano-scale, 
and iv) it evolves in the presence of soft charged inter¬ 
faces. From the point of view of the function, these mate¬ 
rials are also alike. In one case we must understand how 
particular species reach specific active sites on proteins 
or lipidic membranes, and how this relates to the sur¬ 
rounding water. In the other, we aim at grasping the 
mechanisms behind transport of ions within the ionic 
channels, which depend on the shuttling efficiency of the 
water molecules. 

Here, we try to elucidate these issues by investigat¬ 
ing the transport properties of water molecules confined 
in self-organised ionic surfactants phases, across a large 
range of hydration levels. We have selected a perfluo- 
rosulfonic surfactant molecule to insure perfectly phase- 
separated morphologies with well-defined sharp inter¬ 
faces, arising from the extremely hydrophobic character 
of the tail and important super-acidity of the head group. 
We have considered the numerical model for ionic surfac¬ 
tants of Ref. El. which is able to stabilise well-ordered 
self-assembled phases, ranging from lamellae to reverse 
micelles at low and high water content, respectively. We 
have systematically followed the constrained evolution 
of water adsorbed in these mobile soft charged environ¬ 
ments, on extended time scales. On the basis of data 
sets of unprecedented statistical accuracy, we demon¬ 
strate that it is possible to rationalise in a comprehensive 
picture both the average sub-diffusive character of the dy¬ 
namics, due to the subtle interplay between confinement 
at the nano-scale and direct interaction with interfaces 
with different topological properties, and the spatially- 
resolved heterogeneous behaviour in different regions of 


the ionic channels [18]. We also believe that our work 
both highlights the similarities in transport properties of 
systems of biological interest and synthetic materials, and 
provides physical insight for a convincing interpretation 
of phenomenological modelling of experimental data. 


II. THE MODEL 

In Ref. [TT] , we described a coarse-grained numeri¬ 
cal model for (explicit) water confined in ionic surfac¬ 
tants phases. We considered a three-components mix¬ 
ture, formed by the surfactant macro-molecules, the po¬ 
lar solvent (water), and the solvated counter-ions (hydro¬ 
nium). We chose a united-atoms representation for the 
surfactant molecule, inspired by the model for the side- 
chain of Nation of Ref. m- The hydrophobic part of 
the molecule is represented by a chain of seven Lennard- 
Jones (LJ) beads of neutral electric charge, each one 
coarse-graining an entire CF 2 group. Similarly, the head 
group is schematised by two charged LJ beads, one for 
the sulphur atom S and one for the O 3 group, with a 
total charge q = —e. The covalent bonds are modeled by 
quadratic potentials. The non-bonded interaction forces 
were accurately parametrised, in order to generate an ex¬ 
tremely neat phase separation between the hydrophobic 
matrix and the ionic channels, at all investigated water 
contents. The popular SPC/E model was used for the 
water molecules. All details can be found in E3- 

The water content A, i.e., the number of water 
molecules per ionic surfactant macromolecule, influences 
the macroscopic self-assembly and generates a sequence 
of phase changes, encompassing lamellar, cylindrical and 
micellar phases by increasing hydration. In Fig. [I] we 
show typical snapshots of the simulated systems, at 
the indicated values of A, increasing from left to right. 
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Also, in Fig. [2] we show the A-variation of the aver¬ 
age confining size and the surface of the generated hy¬ 
drophilic/hydrophobic interfaces. (This last quantity 
has been calculated by determining a polyhedral mesh 
around the beads pertaining to the surfactants |2Dj, and 
is normalised to the total volume it contains.) Both 
data sets generally increase with A following, however, 
non-trivial patterns which also highlight the boundaries 
between the different phases. The behaviour of the 
confining size can be attributed to the system swelling 
upon increasing water loading. A more subtle argument, 
based on an increasing curvature radius of the interfaces, 
is needed to fully explain the details of the interfaces 
growth PH- 

The interest of this model in the present context is 
two-fold. First, it constitutes a very effective tool to 
autonomously generate very different confining environ¬ 
ments, modulating both the typical confinement sizes 
and the extent and topological features of the formed 
interfaces. Second, the extremely neat phase separation 
provides us with the possibility of precisely identify and 
characterise the hydration layer adjacent to the interface 
with the hydrophobic matrix, complementary to the vol¬ 
ume occupied by water of relatively higher bulk-like char¬ 
acter, far from the boundaries. Both points are crucial 
for rationalising the dynamical behaviour of the adsorbed 
fluid, as it will be clear below. 
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FIG. 2. Variation with the water content, A, of the av¬ 
erage size of the self-assembled ionic channels and the ex¬ 
tension of the generated hydrophobic/hydrophilic interfaces, 
normalised to the volume occupied by the surfactant macro¬ 
molecules. Data for the confining size are similar to those 
reported in El, re-evaluated from the new system snapshots. 
The arrows indicate the axis relevant for each data set. 


IV. RESULTS 


III. SIMULATION DETAILS 


We have used the parallel high-performance computing 
tool for Molecular Dynamics simulations LAMMPS [2Tj 
to generate the trajectories. In order to assure an opti¬ 
mal thermalisation of the systems at all hydrations, we 
reinitialised the simulation boxes starting from the self- 
assembled configurations analysed in HZ], We consid¬ 
ered 9 values of the hydration A ranging in the inter¬ 
val [1,32], which encompass all interesting regions of the 
phase diagram. A pure water system (bulk) trajectory 
was also generated, as a reference. Note that the mobil¬ 
ity of the surfactant macro-molecules lifts the need for an 
additional average over multiple realisations of the con¬ 
fining matrix. We verified in a few cases, however, that 
equivalent results are found when preparing the system 
in different initial conditions. 

At each A, we produced a trajectory spanning a to¬ 
tal time r max = 10 ns, a dynamical window sufficiently 
extended to allow comfortably observing the very slow 
dynamics in the most confined cases. The time step used 
for the numerical integration of the equations of motion 
is St = 1 fs. We worked in the constant temperature 
and pressure (NPT) ensemble, with T = 300 K and 
P = 1 atm, with the characteristic time of both ther¬ 
mostat and barostat set to tt = Tp = 1 ps. We dumped 
complete systems configuration with a frequency of 10 ps, 
for subsequent analysis. 


In this Section we describe our computer simulation 
data. We have fully characterised the dynamics of wa¬ 
ter molecules adsorbed in the surfactants nano-phases 
at the different hydration levels, in terms of mean- 
squared displacements and intermediate scattering func¬ 
tions. We have next completed our analysis by probing 
well-designed correlation functions associated to the typ¬ 
ical life-times of water molecules in the hydration layer, 
in contact with the confining hydrophobic matrix. 


A. Evidence of sub-diffusive behaviour 

By using the produced system trajectories, we first 
computed the water molecules mean-squared displace¬ 
ment (MSD) 
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where N w is the number of water molecules, and r.;(f) 
is the position vector of the oxygen atom pertaining to 
molecule i at time t. In Fig. [3] (top) we show our data 
at the investigated various hydration levels (higher on 
top), together with the results for bulk water. We ob¬ 
serve two interesting features: first, by decreasing the 
hydration, the mean-squared displacement at fixed time 
decreases, as generically expected due to an increasing 
crowding associated to the surfactant units. Second, and 
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FIG. 3. Main panel: Mean-squared displacements (r 2 (f)}, 
calculated from the time-dependent coordinates of the wa¬ 
ter oxygen atoms. The solid lines are the power-law fits 
described in the main text, the dashed lines indicates the 
limiting Fickian (top) and highly sub-diffusive (bottom) lim¬ 
its, corresponding to a = 1 and 0.4, respectively. Inset: X- 
dependence of the exponent a, determined from (r 2 (t)) oc t a . 


more intriguing, while the data-points for the bulk fall 
on a straight line of unitary slope (the linear depen¬ 
dence expected for water at ambient conditions), this 
slope progressively decreases with the hydration, indi¬ 
cating a power-law behaviour with an exponent less than 
one. This is characteristic of sub-diffusion [52 . 

To quantify this feature, we fit the MSD at long 
times to a power-law of the form ( r 2 (t )) oc t a , with 
0 < a = a(A) < 1. We show the results of the fit as 
lines in Fig. [3j and the values of the exponent a(A) in 
the inset of the same figure. We observe that a varies 
slowly from a value close to 1 at A = 32, correspond¬ 
ing to the expected Fickian diffusion, to a ~ 0.85 at 
A = 5, close to the boundary of the lamellar phase. This 
behaviour seems to correlate with the corresponding de¬ 
crease of the available interface shown in Fig. [2] Next, 
a abruptly decreases to a value ~ 0.4 at A = 2, cer¬ 
tainly controlled by the extreme confining sizes at al¬ 
most constant available interface. Note that this value 
is even lower than a = 0.5, corresponding to the case 
of single-file diffusion na, where molecules are organised 
in uni-dimensional configurations, with a vanishing posi¬ 
tion swap probability for adjacent molecules. Altogether, 
these data demonstrate that our model allows to contin¬ 
uously tune the degree of sub-diffusion by varying the 
hydration level, interpolating between the extreme cases 
of Fickian and single-file-like diffusion. 

In this context, the fluctuation-dissipation Einstein re¬ 
lation relating the mean-squared displacement to the dif¬ 
fusion coefficient, does not hold in the usual form. We 
can, however, define a time and length-scale dependent 


generalised (effective) diffusion coefficient, D e g , from the 
local slope of the mean squared displacement as [24], 

d t (r 2 (t)) = j(r 2 (t))=6D eS . (2) 

The effective diffusion coefficients and the corresponding 
( r 2 (t)) are plotted parametrically in Fig. [I] (top). Con¬ 
sistently with the data of Fig. [3] the most diluted cases 
exhibit diffusion coefficients that are rather constant and 
close to the value of the bulk. In contrast, upon decreas¬ 
ing hydration, U e // is strongly suppressed, with this ef- 
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FIG. 4. Top: The effective diffusion coefficient, D a g , com¬ 
puted from Eq. ([ 2 ]), plotted parametrically against ( r 2 (t )). 
The solid black line is the result in the bulk for the considered 
water model, corresponding to D ~ 2.74 x 10” 5 cm 2 /s. The 
dashed vertical lines, indicate the (squared) length scales ~ 6 
and ~ 1000 A 2 used to select the values of D e g shown in the 
bottom panel, as discussed in the text. Bottom: The values 
for D e a selected by the two verticals lines of the above panel 
are shown as closed symbols. We also show the diffusion coef¬ 
ficients obtained by PFG-NMR and QENS spectroscopy [231 . 
which probe the dynamics on long and short length scales, 
respectively. 







































5 


feet becoming increasingly evident as the measured (r 2 ) 
increases. 

An important remark is in order at this point. The 
data shown in Fig. [4] (top) provide us with a unique ratio¬ 
nalisation of the diffusion properties probed by different 
experimental spectroscopic techniques. We illustrate this 
point in Fig. [4] (bottom), where we compare our results 
with the values obtained by Pulsed-Field Gradient Nu¬ 
clear Magnetic Resonance (PFG-NMR ) [23’ and Quasi- 
Elastic Neutron Scattering (QENS) experiments on gels 
of PFOSA surfactants [23]. We recall that PFG-NMR 
(open circles) measures the self-diffusion of water on large 
distances, of the order of the micrometer (/ im ), with 
a typical time scale of a few milliseconds, while QENS 
(open squares) probes proton displacements of the order 
of a few angstroms (A) , on time scales ranging from 1 to 
typically a few hundreds picoseconds. 

Following an uncommon protocol, we sample our data 
at fixed values of (r 2 ) = (l*) 2 , compatible with these 
length scales (vertical dashed lines in the top panel of 
Fig. [4]). Obviously, the /im range explored by PFG-NMR 
is far beyond the limit of our simulations. We note how¬ 
ever, that at least for the highest hydration levels, D e g 
is almost length-scale independent for ( r 2 (t)) > 10 3 A 2 . 
This allows us to safely fix an observation threshold at 
l* ~ 32 A. Similarly, we have determined a short length- 
scale value of D e s with l* ~ 2.4 A, which is compatible 
with the typical nano-scale distances probed by QENS. 
Interestingly, the simulation points nicely reproduce the 
hydration dependence of both experimental data sets, 
without resorting to any additional consideration about 
different dynamical processes on different time scales. 
Our model therefore seems to provide not only an ade¬ 
quate dynamical behaviour compared to the real system, 
but also a framework to interpret experimental results, 
without needing any fitting procedure (beside the deter¬ 
mination of a). In the following we will demonstrate 
that it also is an efficient tool to pinpoint the origin of 
the sub-diffusive behaviour described here. 


B. Existence of two dynamical populations 

Additional information on the dynamics of the water 
molecules is carried by the (self-)intermediate scattering 
function, 


F s {q,t) 


1 N„ 

_L V />«J-(ri(t)-ri(0))\ 

n w ^\ 4w 


( 3 ) 


where ()|g|= 9 is the spherical average over wave vectors 
of modulus q. A strong connection exists in the bulk 
between the F s (q,t) and the diffusion coefficient. We 
expect it to partially survive, even if in some modified 
form, in confinement conditions. The Fourier transform 
of Eq. i> is the dynamical structure factor, S(q, w), mea¬ 
sured in the QENS experiments mentioned above. Here 
we focus on the value q* = 0.7 A -1 [26| . We show our 
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FIG. 5. The self-intermediate scattering functions, F s (q,t), 
calculated from the oxygen atoms of the water molecules via 
Eq. ([ 3 ]) (symbols). These data correspond to q* = 0.7 A -1 , at 
the indicated values of A. We also show as solid lines the best- 
fits to the data, according to Eq. Q. These results strongly 
support the hypothesis of the existence of two dynamical pop¬ 
ulations, as discussed in details in the main text. 


data at all the indicated values of A in Fig. [5] together 
with the reference case of the bulk. We observe that the 
decay of F s (q,t) becomes progressively slower as the hy¬ 
dration is decreased, and that the general shape of these 
functions is quite complex. In the following we quantify 
these observations. 

The long-time behaviour of F s (q,t) for bulk molecular 
liquids is often represented as a stretched exponential of 
the form oc . Here r is the typical time scale as¬ 

sociated to the observed relaxation, and 0 < (3 < 1 is the 
stretching parameter. The case /3 = 1 corresponds to the 
simple exponential behaviour, /3 < 1 describes the case 
where non-trivial processes concur to the average relax¬ 
ation of the correlation function. As one could expect 
from the extremely complex behaviour of (r 2 (f)), this 
simple functional form alone is not able to accurately re¬ 
produce our data, especially at low A where the degree 
of sub-diffusivity is higher. 

The simplest alternative choice is to assume the exis¬ 
tence of two molecular populations characterised by dif¬ 
ferent dynamical properties [23] [27[ [28]. In this quite 
crude representation, a subset of the total number of wa¬ 
ter molecules evolve diffusively, not very differently than 
those in bulk conditions, on a relatively fast time scale. 
The remaining part, in contrast, include all the complex 
effects coming from confinement in the surfactant phase, 
evolving on slower time scales. The precise features of 
the two populations should intuitively depend on the de¬ 
tails of the confinement environment and, ultimately, on 
A. Given the additivity of Eq. ([3]) , we can therefore con¬ 
sider the very simple model, 

Fs(q, t) = /sd e ^ t sd) + (1 -/ SD ) e (m), 


( 4 ) 
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where a simple and a stretched exponentials encode the 
dynamics of the diffusive (D) and sub-diffusive (SD) wa¬ 
ter molecules, respectively [ 21 ] ■ Since this latter popu¬ 
lation is the slowest, we can expect td < ( tsd ), where 
the average time scale associated with the sub-diffusive 
population is (tsd) = {tsd/P) T (1//3), with r(x) the 
Gamma function. We can also obtain the fraction of dif¬ 
fusive molecules, as fn = (1 — Jsd)- 

We were able to obtain excellent fits to our data with 
this model, that we show as solid lines in Fig. [5] The 
obtained values of fsD span the entire possible variation 
range, assuming a value of 0.05 at the highest A and of 
0.95 at the lowest available hydration level. This varia¬ 
tion is followed by /3-values ranging from 0.6 to 0.3. These 
observations point to an increasing modification of the 
dynamics compared to the bulk unconstrained case, and 
to a very probable connection between the sub-diffusive 
average character of the dynamics and the existence of 
molecular populations, characterised by heterogeneous 
dynamical properties. We postpone a much more elab¬ 
orated discussion of the calculated values for both the 
time scales and fsD to the following (see Fig. [9]). 

C. The interface as the source of sub-diffusion 

The self-assembled confining environments considered 
here are characterised not only by variable size of the 
ionic domains available for transport (Fig. [2]), but also by 
the variable topological features of the formed interfaces, 
which evolve from the almost perfect planar configura¬ 
tions of the lamellar phase, to configurations of increasing 
convexity in the cylindrical and micellar phases. It is nat¬ 
ural to expect that beside the obvious role played by the 
extremely tight confining sizes at very low A, these well- 
developed interfaces significantly perturb the dynamics 
of the adsorbed fluid molecules, participating to the gen¬ 
eration of the observed sub-diffusion. More specifically, 
the surfactant sulfonic acid groups are hydrophilic and 
negatively charged, resulting in a strong attractive po¬ 
tential between the interface and fluid molecules within 
the range of these interactions. This leads to the hy¬ 
pothesis that the fluid molecules in the vicinity of the 
interface can be trapped at the interface, which reduces 
their possibility of diffusing at the normal rate. This 
constrain would next be transferred to successive fluid 
layers, with a strength reduced as the distance from the 
interface increases. 

On this basis, we now focus on the dynamical prop¬ 
erties of the fluid adsorbed at the interface only , by se¬ 
lecting the sub-set of fluid molecules comprised in the 
first hydration layer of the surfactant heads. We there¬ 
fore define the binary presence function, p{ r), for a water 
molecule at position r: p(r) = 1 if it exists a charged 
surfactant head at r 0 , such that |r 0 — r| < r c , and zero 
otherwise. In what follows we have fixed r c = 4 A. Ap¬ 
plying this operator to the position of the oxygen atoms 
for every snapshot of our trajectory, we have obtained N w 
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residence time of H 2 O at the interface [ps] 
1+ 4 x 8 * 12 □ 32 ■ 



FIG. 6. Top: Probability distributions of the residence times 
observed during the 10ns simulation run. The long-time be¬ 
haviour is reminiscent of a power-law, and in the inset we plot 
the exponent of power-law fits. Bottom: Distribution of the 
residence times of the water molecules present at the inter¬ 
face a the end of the simulation run. The black arrow points 
to the long-time bump in the distribution that is particularly 
visible at low hydration (A < 6). 

series of N s binary values, {pi, n }- Here, N s is the number 
of system snapshots, each corresponding to a simulation 
time t n = n St, with n = 0,..., N a — 1, and i = 1,..., N w . 

From the lists {pi, n } we can in principle determine the 
probability distribution of the residence times at the in¬ 
terface of all molecules, by simply counting the number of 
non-zero segments of given size (for each i). This should 
allow us to define an average survival time at the inter¬ 
face. Unfortunately this is problematic, as it is clear from 
the data sets shown in Fig. [ 6 ]top, at the indicated values 
of A. These distributions show significant power-law tails 
oc t~ x (solid lines in the main panel and points in the in¬ 
set), with even x < 2 at the lowest hydrations (A < 2) 
which hinders a consistent definition of average values, 
and 2<x<3for2<A<10, which results in an un- 
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defined variance. This also implies that there is always 
a significant probability to observe life-times correspond¬ 
ing to the total time length of the simulation, making 
extremely hard adequate statistical sampling. This con¬ 
clusion is verified from the data of Fig. [6] (bottom), where 
we plot the same probability distributions as above, with 
the additional constraint that the considered non-zero 
segments must include r max (3(1;. Clearly, at low hy¬ 
dration values the most likely residence time of water 
molecules present at the interface at the end of the run 
is the duration of the run, as indicated by the arrow. 

These observations are of particular interest, because 
long-tailed distributions of survival times are known to 
give rise to sub-diffusive random walks. This is the 
case, for instance, of the continuous-time random walk 
(CTRW), where particles are allowed to move after hav¬ 
ing been arrested for a waiting time, t w , drawn from a 
power-law distribution. Evidences of CTRW have been 
observed in simulations of water confined in lipid bi¬ 
layers mm, among other cases. 

The difficulties described above can be partially cir¬ 
cumvented by determining the time scale associated to 
the variation of well-designed correlation functions. We 
first consider the correlation function of the survival 
times of water molecules laying within the first hydra¬ 
tion layer of the sulfonated groups as 


Ci(t„) 
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/ v -- ' N s — 1—n 
\2^s =0 


(1 Pi,n+s) 0 


n+s— 

m—s 


1 



1 - n/N s 




Here, ()j is the average over all water molecules i = 
1 ,...,N W , and t n = t n+s — t s . Note that the product 
at the numerator is non-zero only if particle i survives 
at the interface at all times t s < t < t n+Sl subsequently 
detaching exactly at t n+s . We plot our data in the main 
panel of Fig. [7] as symbols, at the indicated values of A. 

Already by visual inspection it is clear that the typi¬ 
cal lifetime of water molecules at the interface increases 
by lowering A. Interestingly, we also find that two dif¬ 
ferent mechanisms seem to concur at the relaxation on 
different time scales. Our data can be represented at all 
hydrations by a modified exponential, C\{t) oc 
at short times, crossing-over to a power-law behaviour, 
C'i(t) a at longer times. Note that the cross-over 

point depends on A. The fitting curves are shown as solid 
and dashed lines in Fig. © 

In our understanding, the faster process is associated 
to the exchange mechanism of water molecules which de¬ 
tach from the interface and transit in higher-order hydra¬ 
tion layers, on a (well-defined) typical time scale t s ( A). 
Indeed, the functional form used is suggested by an anal¬ 
ogy with the decay of a subset (population) of a given 
system, characterised by some particular feature. Simi¬ 
larly to the case of radioactive decay, for instance, if the 
probability to disintegrate (or to detach from the inter¬ 
face in our case) is uniform, then C\(t ) ss e ~ t ^ Ts . In con¬ 
trast, if the probability to detach is not a constant, then 



FIG. 7. Main panel: Survival-time correlation functions of 
water at the interface, Ci(t), computed via Eq. (©, at the 
indicated values of A. Insets: A-dependence of the two fitting 
parameters and a a , obtained by representing the data with 
a modified exponential decay at short-times and a long-time 
power-law. The rationale for the fits shown as solid lines in 
the main panel is discussed in details in the main text. 


Ci ( t ) can relax both faster or slower than a simple expo¬ 
nential. The mildest deviation from the pure exponential 
is Ci(f) ~ e ~( t / T s) f>s , with the exponent /3 S quantifying 
the degree of the deviation: for /3 S < 1 (stretched) the 
decay is slower than exponential, faster for /3 S > 1 (com¬ 
pressed). In more extreme cases, scale free power-law 
behaviour have also been observed mm- 

Our analysis suggests that t s increases from 60 to 
22000 ps as the hydration A decreases from 32 to 1 (see 
Fig-©, meaning that the exchange mechanism at the in¬ 
terface is activated by hydration. Also, we have found 
that /3 S varies from 0.6 at A = 1 to 0.8 at A = 32 
(inset of Fig. ©. These observations therefore seem to 
demonstrate that hydration controls both typical time- 
scale and statistical features of the postulated exchange 
mechanism. In addition, the self-similar behaviour ob¬ 
served at longer times clearly is the echo of the power-law 
tails of the static distribution functions of the {pi, n } of 
Fig.© It is natural to expect that this relaxation stems 
from trapping at the interface of water molecules on ex¬ 
tremely extended time scales. This would also directly 
allow to identify interfaces as the main source of water 
sub-diffusion, rationalising this mechanism in terms of a 
CTRW approach. 

It is not possible, however, to associate a well-defined 
time scale to this dynamics in terms of Eq. ©. Indeed, 
in order to highlight the dynamics most controlled by the 
interface, beyond the exchange process quantified above, 
one should lift the excessive constraint of water molecules 
uninterruptedly bounded to the interface during the time 
window t n implicit in that definition. We therefore define 













V. DISCUSSION 



FIG. 8. Time autocorrelation functions of {pi,n} (symbols), 
computed via Eq. , at the indicated values of A. The con¬ 
tinuous lines are the results of the fits to our data according 
to Eq. 0. A comprehensive discussion of these results is 
provided in the main text. 


the time autocorrelation function of {pj,„}, 


c 2 (t n ) 


/^N B -l-n \ 

\Z-^s =0 Pi,s+nPi,sJ . 



( 6 ) 


At variance with the case of Eq. ([5]), with this defini¬ 
tion water molecules that are at the interface at times 
t s and t a +n contribute to C 2 (t n ), notwithstanding their 
faith between t s and t s + n , and at t s+n+ 1 . This func¬ 
tion therefore still selects molecules that detach from the 
interface, but also those that return to the interface, at 
some later time after leaving. We expect, in general, that 
this dynamics takes place on time scales comparable to 
or longer than t s . We plot our data as symbols in Fig. [8j 
at the indicated values of A. 

The data show a two-relaxation pattern, that we de¬ 
scribe as, 


C 2 (t) = ce~^ Pi + (1 - c)e~ { - t/T ^°°. (7) 


We consider here two points. First, the initial decay of 
Gb(t) takes place on a quite large range of time scales r,, 
ranging from ~ 100 ps at A = 32 to ~ 1 ns at A = 1. 
This quantifies dynamics at the interface, as we will see 
below. Second, the functions exhibit an abrupt singular 
behaviour as t approaches the total time window of our 
simulation, r max , and this cutoff is more pronounced at 
lower hydrations. We indicate this limit time scale as 
Too, which characterises molecules that are permanently 
attached to the interface on the time scale of our simu¬ 
lation. We postpone a comprehensive discussion of these 
data to the next Section. 


We now describe the implications of the data reported 
in the previous Section, and provide a consistent unique 
picture of our observations. From the analysis above, 
we have been able to extract five time-scales relevant for 
the confined dynamics of the water molecules. In par¬ 
ticular, T£> and tsd are associated to the relaxation of 
the intermediate scattering function, t s characterises the 
survival-time correlation function C'i(t), r, and t, oo have 
been determined from the binary presence autocorrela¬ 
tion function C 2 (t) . We plot all data as a function of the 
hydration level A in the top panel of Fig. [9] 

We see that the long-time component of the autocorre¬ 
lation of the binary presence, Too, is A-independent, and 
of the order of the total simulation time. We associate 
this time-scale to the water molecules that stay in the 
interfacial layer for times that are compatible with the 
duration of the entire simulation run. As already antic¬ 
ipated, this is an artificial feature which originates from 
arbitrarily fixing a bound from above on the available ob¬ 
servation time. As a consequence, all dynamics on time 
scales larger than Too cannot be probed in equilibrium in 
our simulation. 

Interestingly, the time-scales Tj of the short-time com¬ 
ponent of C 2 (t) and t s of the initial decay of C\ ( t ) are, at 
variance, strongly reduced as the water content increases. 
We can associate these relaxations to an exchange mech¬ 
anism between the interfacial region and the rest of the 
aqueous domains, which is activated by the increase of 
the hydration. This idea is strongly supported by a sur¬ 
prising observation: Tj and t s are very similar, especially 
at low-A, to the time scale tsd that we extracted from 
the F s (q*,t), and associated with the sub-diffusive frac¬ 
tion of water molecules. These therefore seem to be those 
pertaining to the population with a dynamics of primary 
interface-dominated character. 

An other observation additionally strengthen our pic¬ 
ture. In Fig. [ 9 ] (bottom) we show two quantities: the 
static fraction over the total number of water molecules 
lying at the interface, /), which we obtain by simply 
counting the number of molecules comprised in the hy¬ 
dration layer; and the dynamic fraction of sub-diffusive 
particles, fsD, obtained from Eq. Q. The two data 
sets exhibit very similar variation on the whole hydra¬ 
tion range, and the two curves actually coincide at the 
lowest hydration ranges, A < 2. At higher hydration lev¬ 
els, fsD is close to, although always higher than, _/). We 
infer from the similarities of the typical time scales and 
associated fractions (over the total number) that inter¬ 
facial and sub-diffusive populations coincide, or at least 
that the sub-diffusive population is a superset of the in¬ 
terfacial one. 

We propose a schematic overview of our discussion in 
Fig. |10[ where we partition into three main areas the 
surfactant phase environment, with both the confining 
matrix and the complementary space available for trans¬ 
port in the ionic channel, of typical size l: i) the solid- 
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FIG. 9. Top: A-dependence of the relevant time-scales asso¬ 
ciated to the dynamics of the confined water molecules: rn 
and rsD are associated to the relaxation of the intermediate 
scattering function, r s characterises the survival-time corre¬ 
lation function, Ci(f), n and Too have been determined from 
the binary presence autocorrelation function, C 2 (t). A com¬ 
prehensive discussion of these data is included in the main 
text. Bottom: Static fraction of water molecules residing at 
the interface, /;, together with the dynamic fraction of sub- 
diffusive molecules, fsD- The two data-sets coincides at low 
hydration. This point, together with the discrepancy at high- 
A are described in details in the main text. 


like interface (bottom), which comprises the surfactant 
hydrophilic heads; ii) the hydration layer (centre), that 
covers the interface and comprises all water molecules 
within r c = 4 A from the the nearest surfactant head; 
and in) the bulk region (top) consisting of all remaining 
water molecules and extending over a typical length scale 
lb, which depends on A. Obviously, l = 2(4 + 4/2). 

The hydration layer is characterised by the fraction of 
molecules it comprises, /,, together with the characteris¬ 
tic time of their interface-dominated dynamics, r,. The 
fraction of water at the interface coincides at low-A (and 
has a very large overlap at high-A) with the fraction of 


sub-diffusive water (/,; ss Jsd), and one finds r, ss tsd- 
The molecules at the boundary with the region extending 
beyond the first hydration shell and comprising the re¬ 
maining fraction, 1 — /,; ss _/4>, participate to the exchange 
process with the molecules pertaining to the bulk-like 
reservoir. This takes place on a time scale r s , which also 
correlates to the sub-diffusive time scale of the interme¬ 
diate scattering function, with t s ~ tsd- 

The sub-diffusive character of the mean-squared dis¬ 
placements (Fig. [3]) also finds a clear rationalisation 
in this picture. (r 2 (t)) grasps an average dynamics, 
which results from a A-dependent interplay of interface- 
dominated and bulk-like behaviour. The source of the 
observed sub-diffusion is trapping of water molecules at 
the interface, with a power-law distribution of trapping 
times, as demonstrated in Fig. [7] This breaks the va¬ 
lidity of the central limit theorem in a fashion remi¬ 
niscent of a continuous time random walk, immediately 
providing a rigorous justification of the observed non¬ 
linearity of ( r 2 (t )) at long times. Also note that in gen¬ 
eral fi < fsD, and at high hydrations (above A = 6) 
t s and Ti are slightly larger than tsd- This indicates 
that water molecules exiting the hydration shell retain 
a certain degree of sub-diffusivity, with this anomalous 
character progressively disappearing when moving away 
from the interface. 


bulk-like region: 1-fj 
diffusive water: f D ~l-f i/ t d 


=5— 

Ts Tsd hydration layer: f^ 
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sub-diffusive water, f ND —f,, TND'""' T i~ T s 
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solid-liquid interface 


FIG. 10. Schematic overview of the interplay between the 
structure of the surfactant phase and the transport properties 
of the confined water molecules. This sketch complements the 
detailed discussion contained in the main text. 


Confinement size, L=2(L h +L b /2) 
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FIG. 11. Joint probability distribution of the total displacement, l, of water molecules on the entire time scale of the run, r m ax, 
and the total time, t h , they have spent (even discontinuously) in the hydration layer of a surfactant head, on the same time 
window. The probability is represented as colour maps, for the indicated values of A: yellow indicates high probability, indigo 
lower. In the insets we show the projection on the xy -plane of the trajectory followed by two molecules chosen at random at 
the corresponding A, with the incremental distance travelled at each time indicated by the displayed colour code. 


To further clarify this point, we have quantified 
the correlation between the total displacement, lj = 
Y^n^i 1 l r i(^n) — r i (t n ~ 1)1 > of water molecule i on the 
entire time scale of the run (r max ~ 10 ns) and the to¬ 
tal time, t%, it has spent (even discontinuously) in the 
hydration layer of a surfactant head, on the same time 
window. In Fig. El we show the joint probability distri¬ 
butions of these two quantities as colour maps, for the 
indicated values of A; yellow indicates high probability, 
indigo lower. We observe a clear trend as the hydration 
increases. 

At very low water content (A = 1), almost all water 
molecules remain in the hydration layer for most of the 
duration of the simulation run (th ~ T max ), and the total 
distance travelled does not strongly depend on the time 
spent at the interface, fluctuating around l ~ 750 A. In 
the inset we show the projection on the xy -plane of the 
trajectory followed by two molecules chosen at random, 


with the distance from the nearest O 3 bead indicated by 
the displayed colour code. The motion is clearly very con¬ 
strained (probably bi-dimensional) in both cases, with a 
very limited portion of space explored. For A = 4 and 8 , 
however, the probability to observe very high values of 
th progressively decreases, while the distribution of the 
th becomes wider. This effect is accompanied by longer 
l travelled at shorter th, with a degree of variation in¬ 
creasing with A, as confirmed by the data in the insets. 
At A = 32, eventually, the probability accumulates in 
the high-//low-f/j region, indicating that practically all 
molecules transit between the interface and the bulk-like 
region (the trajectory of one molecule is shown in the 
inset, and clearly shows segments of high-speed travel 
far from the interface, and “sticky” regions that trap the 
molecule). 

The conclusion arising from the above arguments is 
that the picture based on two well-identified dynami- 
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cal population, that fully accounts for the behaviour in 
the low-A regime, must be modified at high water con¬ 
tent, where we observe an increasing rate of the exchange 
mechanism and of the number of molecules participating 
to it. Indeed, it is not possible to account for the complex 
distributions of Fig.[ll]in terms of an average r, only, and 
one should rather consider a A-dependent set of popula¬ 
tions, characterised by a distribution of t % ( tnd )• 

VI. CONCLUSIONS 

In this work we have followed numerically the dy¬ 
namics of water molecules constrained in self-organised 
ionic surfactant phases. By varying the water content, 
we have considered an extended range of soft mobile 
confining environments, characterised by variable con¬ 
fining scales and topological properties of the generated 
hydrophilic/hydrophobic interfaces. These two features 
contribute in a highly non-trivial fashion to transport of 
the adsorbed fluid, that we have managed to characterise 
in terms of a space-dependent dynamical behaviour. 

More in details, we have demonstrated that both the 
observed sub-diffusive character of the average mean 
squared displacements, increasing by decreasing the 


length scale of the confinement, and the details of the 
relaxation of the intermediate scattering function, can 
be rationalised in terms of two dynamical populations. 
The first is predominantly interface-dominated, with a 
marked non-diffusive behaviour, the second has a bulk¬ 
like character. The average dynamical behaviour of the 
adsorbed fluid is primarily determined by the exchange 
mechanism between the two dynamical groups, whose 
precise features depend on the water content. This sim¬ 
ple picture allows for a complete quantitative character¬ 
isation of dynamics at low-A, while at high-A we have 
proposed that an entire A-dependent set of populations 
must be considered to obtain the same accuracy. 

Overall, we have managed to provide, to the best of 
our knowledge, the first comprehensive picture of the dy¬ 
namics of water molecules in soft confining environments, 
ranging from the case of extremely tight channels, to that 
where the interface is in contact with an infinite reser¬ 
voir of bulk water. These two limits cover, we believe, 
most part of the typical situations encountered in hy¬ 
drated natural and synthesised materials. Our model, 
as a consequence, could play for soft materials a critical 
paradigmatic role similar to that of carbon nanotubes for 
confinement in rigid ordered matrices. 
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